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ABSTRACT 

A two-dimensional cloud-resolving model is used to examine the possible role of gravity waves generated 
by a simulated tropical squall line in forcing the quasi-biennial oscillation (QBO) of the zonal winds in the 
equatorial stratosphere. A simulation with constant background stratospheric winds is compared to simulations 
with background winds characteristic of the westerly and easterly QBO phases, respectively. In all three cases 
a broad spectrum of both eastward and westward propagating gravity waves is excited. In the constant background 
wind case the vertical momentum flux is nearly constant with height in the stratosphere, after correction for 
waves leaving the model domain. In the easterly and westerly shear cases, however, westward and eastward 
propagating waves, respectively, are strongly damped as they approach their critical levels, owing to the strongly 
scale-dependent vertical diffusion in the model. The profiles of zonal forcing induced by this wave damping 
are similar to profiles given by critical level absorption, but displaced slightly dow nward. The magnitude of the 
zonal forcing is of order 5 m s 1 day '. It is estimated that if 2% of the area of the Tropics were occupied by 
storms of similar magnitude, mesoscale gravity waves could provide nearly 1/4 of the zonal forcing required 
for the QBO. 


1. Introduction 

In their theory of the quasi-biennial oscillation (QBO) 
of the zonal mean winds in the equatorial stratosphere, 
Lindzen and Holton (1968, hereafter LH68) argued that 
zonal forcing induced by a continuous spectrum of long- 
period vertically propagating gravity waves was re- 
sponsible for the downward propagating easterly and 
westerly mean wind regimes that characterize the QBO. 
In LH68 the wave-induced forcing of the mean flow 
was attributed to absorption of these vertically propa- 
gating waves at "critical levels” where their Doppler- 
shifted horizontal phase speeds vanish. According to 
this hypothesis, waves with a broad spectrum of zonal 
phase speeds would be absorbed in a shear zone in which 
the mean zonal wind speeds match the phase speeds of 
the waves. LH68 showed that if the gravity wave spec- 
trum contained components with both positive and neg- 
ative zonal phase speeds, the critical level absorption 
mechanism could lead to a mean zonal flow oscillation 
with alternating downward propagating westerly and 
easterly shear zones, provided that a mechanism existed 
to initiate the shear zones in the upper stratosphere. In 
LH68 this initiation mechanism was taken to be the 
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observed semiannual oscillation centered near the equa- 
torial stratopause. 

Somewhat later, and motivated by observational stud- 
ies of large-scale equatorial waves (Maruyama 1969; 
Wallace and Kousky 1968), Holton and Lindzen (1972, 
hereafter HL72) offered an alternative to the LH68 mod- 
el in which the zonal forcing was attributed to two large- 
scale monochromatic equatorial wave modes, an east- 
ward propagating Kelvin wave, and a westward prop- 
agating Rossby— gravity wave. In HL72 the wave-in- 
duced forcing of the mean flow was attributed to 
radiative damping of the long-period Kelvin and Ross- 
by-gravity waves rather than to critical level absorption. 
Since the rate of radiative damping increases rapidly as 
the Doppler-shifted phase speed decreases, Kelvin 
waves will be preferentially damped and produce east- 
ward accelerations in westerly shear zones, while Ross- 
by-gravity waves will be preferentially damped and pro- 
duce westward accelerations in easterly shear zones. Be- 
cause waves are damped below their critical levels, this 
mechanism no longer requires a separate process to ini- 
tiate shear zones at high levels (Plumb 1977). 

The HL72 model was able to produce an oscillation 
similar in character to that of the LH68 model. However, 
it was soon recognized by Lindzen and Tsay (1975) that 
an additional equatorial easterly gravity wave mode was 
required to account for the observed easterly accelera- 
tion of the QBO. More recently, Dunkerton (1991) 
showed that when vertical advection by the mean up- 
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welling in the equatorial stratosphere is taken into ac- 
count, observed Kelvin and Rossby-gravity wave fluxes 
are much too small to generate the observed QBO. Stud- 
ies by Gray and Pyle (1989) with a two-dimensional 
model and Takahashi and Boville (1992) with a three- 
dimensional model confirmed that in the presence of 
mean upwelling the Kelvin and Rossby-gravity wave 
fluxes must be at least a factor of 2 greater than are 
observed in order to produce a realistic QBO. 

Recently, Dunkerton (1996) considered the QBO mo- 
mentum problem in some detail. He showed that if a 
suitable spectrum of mesoscale gravity waves is added 
to the observed Kelvin and Rossby-gravity waves, the 
combined zonal forcing can produce a realistic QBO 
(including the observed asymmetries between the east- 
erly and westerly shear zones) even in the presence of 
mean upwelling. Dunkerton used the gravity wave mo- 
mentum flux parameterization of LH68, and a Kelvin 
and Rossby-gravity wave damping parameterization 
similar to that used in HL72 (Lindzen 1971). Thus, it 
now appears that a combination of the models of LH68 
and HL72 is required to explain the observed QBO. 

Although LH68 emphasized the likely importance of 
equatorially trapped internal gravity waves of large hor- 
izontal and small vertical scale, their parameterization 
for the wave-induced force was actually derived from 
consideration of internal gravity waves in the absence 
of rotation and, thus, was more applicable to mesoscale 
gravity waves. At the time of the writing of LH68 there 
was little direct observational evidence of such gravity 
waves in the equatorial atmosphere. However, in recent 
years analysis of high-resolution radiosonde data (e. g., 
Sato et al. 1994; Tsuda et al. 1994; Allen and Vincent 
1995; Karoly et al. 1996), data from high-altitude re- 
search aircraft flights (Alexander and Pfister 1995; Pfis- 
ter et al. 1993), and satellite cloud images (Bergman 
and Salby 1994) all suggest a substantial role for con- 
vectively generated mesoscale gravity waves in the 
equatorial region. 

The morphology of gravity waves generated by spe- 
cific equatorial convective systems is, however, still un- 
clear. Little is known about the spectral distribution of 
waves generated by tropical convection, their propa- 
gation into the stratosphere, and their interaction with 
the mean flow. In this work we employ a two-dimen- 
sional numerical model of a tropical squall line. Al- 
though convective systems often show a tendency to 
occur in lines or bands, the squall line is a less general 
case, chosen because its symmetry makes a two-di- 
mensional modeling approach reasonable and because 
the tropospheric properties of such systems have been 
studied in detail elsewhere (e. g., Houze 1993). 

The object of the present study is to examine the 
gravity wave fluxes in the stratosphere simulated by this 
special type of tropical convection and to show that the 
vertically propagating gravity waves generated by such 
a system may play a significant role in the momentum 
budget of the lower equatorial stratosphere. 


2. Numerical model 

The tropical squall line simulations were carried out 
with the nonhydrostatic cloud resolving model of Dur- 
ran and Klemp ( 1 983), who provided details of the cloud 
microphysical parameterizations and the numerical so- 
lution method used in the model. The same model was 
previously used by Fovell et al. (1992) and Alexander 
et al. (1995) to study gravity waves produced by a sim- 
ulated two-dimensional midlatitude squall line. 

The version of the model used here is only slightly 
altered from that used by Alexander et al. (1995). The 
model domain is 900 km wide and 30 km deep. Grid 
spacings are 1.5 km in the horizontal and 0.25 km in 
the vertical. The high resolution in the vertical is de- 
signed to resolve the shallow cold pool near the surface 
and to resolve short vertical wavelength gravity waves 
in the stratosphere. The radiation boundary condition of 
Durran (1995) is used at the upper boundary. In addi- 
tion, to avoid spurious cooling near the upper boundary, 
it is necessary to include a “sponge” layer in which the 
potential temperature is relaxed toward the initial profile 
by Newtonian cooling at a rate increasing uniformly 
from zero at the 27-km level to a maximum of 5. OX 1 0 5 
s 1 at the top boundary. A radiation boundary condition 
is also employed at the lateral boundaries. But, to pre- 
vent noise from developing near these boundaries in the 
troposphere, it is also necessary to include a second- 
order smoothing operator, which is applied to all fields 
at grid points within 30 km of both lateral boundaries. 

To prevent nonlinear instability, weak fourth-order 
diffusion is applied to all prognostic fields in the hor- 
izontal and to the velocity and pressure fields in the 
vertical. The fourth-order vertical smoothing has an e- 
folding decay time of about 400 s for a vertical wave- 
length of 1 km, 6300 s for a vertical wavelength of 2 
km, and has a negligible effect for waves of vertical 
wavelength greater than 4 km. Because the vertical 
wavelength shrinks rapidly as a wave approaches its 
critical level, the fourth-order vertical smoothing has 
important implications for damping of waves in shear 
zones in the simulations reported here. 

In addition to the fourth-order smoothing required to 
prevent nonlinear instability, there is also a parameter- 
ized, Richardson-number-dependent, subgrid-scale mix- 
ing in the model as described in Durran and Klemp 
(1983). For the cases reported here, subgrid-scale mix- 
ing is only important in the troposphere. The fourth- 
order numerical diffusion generally damps waves in the 
stratosphere before the Richardson number becomes 
sufficiently small to activate the explicit subgrid-scale 
mixing parameterization. 

The initial thermodynamic and tropospheric wind 
soundings (Fig. 1) are specified to approximately sim- 
ulate conditions of the West African squall lines studied 
by Redelsperger and Lafore (1988) and Lafore and Mon- 
crieff (1989). The specified thermodynamic sounding 
has an unusually high conditional instability for the 
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Fig. 1. Initial thermodynamic and wind soundings used in the trop- 
ical squall line simulations. The left panel shows temperature (thick) 
and dewpoint temperature (thin). The right panel shows the tropo- 
spheric wind profile used in all simulations. The winds are a constant 
-16 m s _1 above the tropopause in the control simulation. Strato- 
spheric shear cases are shown in Fig. 2. 


Tropics. Simulations with varying thermodynamics 
soundings (not shown) demonstrate that for the two- 
dimensional version of the cloud model used here, rather 
large values of conditional instability are required for 
development of sustained mesoscale convection. 

In order to analyze the possible contribution of grav- 
ity waves generated by squall lines to acceleration of 
mean winds in the QBO, three simulations are com- 
pared. In the first case (labeled the “control”), winds 
above the tropopause are constant with height at —16 
m s _! . The model domain also translates at a constant 
speed of -16 m s _1 . Thus, stratospheric waves are 
viewed in an intrinsic frame of reference (i.e., phase 
speeds relative to the model domain are the Doppler- 
shifted phase speeds in the stratosphere). Since the speed 
of motion of the squall line is only about —5 m s~ l in 
the translating domain, the squall line remains well with- 
in the domain throughout a 12-h simulation period. The 
other two cases (labeled the “easterly shear” and “west- 
erly shear” cases, respectively) have the wind profiles 
above the tropopause used by Holton (1972) to represent 
the easterly and westerly phases of the QBO, respec- 
tively (Fig. 2). 

The simulations for each of the three cases were run 
for 12 h, and field variables were saved at 3-min inter- 
vals for the last 8 h of each simulation. The overall 
structure of the control case at 6 h is shown in Fig. 3, 
which may be compared with Fig. 2 of Alexander et al. 
(1995). Temporal variability in the strength of the squall 
line is illustrated in Fig. 4, which shows a longitude- 
time cross section of rainfall rate. It is clear from Fig. 
4 that the squall line drifts slowly westward in the model 
domain. But the average speed of drift is only —4.9 m 
s -1 in the moving frame, so that the convection remains 
far from the boundaries throughout the integration pe- 
riod. Figure 4, which may be compared with Fig. 3 of 
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Fig. 2. Stratospheric wind profiles for simulations of the easterly 
(dashed line) and westerly (dotted line) phases of the QBO. 

Fovell et al. ( 1 992), indicates that the variability in these 
tropical simulations is less regular than the quasiperiodic 
variability seen in the midlatitude cases previously sim- 
ulated. Thus a longer period was deemed necessary for 
diagnosis of the wave spectra. 

3. Wave momentum fluxes 

The properties of gravity waves produced in the sim- 
ulations can be conveniently characterized using the 
two-dimensional (x, t) spectral analysis approach de- 
scribed in detail in Alexander et al. ( 1 995). This analysis 
is carried out for the period 4 h < ( < 11 h. Because 
average travel times for gravity wave energy across the 
vertical stratospheric domain are 2-3 h, this analysis 
period excludes start-up artifacts. Also, transience will 
be effectively averaged out and will not artificially ap- 
pear as vertical variations in the momentum flux. Only 
grid points in the range 0 km < x ^ 768 km are included 
in the spectral analysis in order to give approximately 
equal weighting to the regions upstream and down- 
stream of the squall line gust front. 

Spectra of vertical velocity and of momentum flux 
show similar features as those discussed in Alexander 
et al. (1995) for a midlatitude case. In that previous 
work, spectral power was concentrated along a line of 
constant vertical wavelength of 8 km, corresponding to 
the depth of the convective heating in the troposphere. 
In this tropical case, there is more variability in the depth 
of convective heating as a function of time. The spec- 
trum of stratospheric gravity waves correspondingly 
shows power scattered more broadly between vertical 
wavelengths of 2-15 km — the upper limit of this range 
associated with the depth of the troposphere. Two pro- 
nounced peaks at 5 and 6 km correspond to the heating 
associated with the rearward propagating cells, which 
have 5-6 km depth in this simulation. Unlike the mid- 
latitude case, no pronounced peaks appear in the fre- 
quency spectrum, suggesting a less important role for 
the “mechanical pumping” mechanism of wave gen- 
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Fig. 3. The control simulation at 6 h of simulation time. Shading represents contours of 
vertical velocity. (Contrast has been enhanced to show the qualitative structure; the full range 
of vertical velocities is ±5 m s '.) Thin lines are isentropes (at 10-K intervals), and the thick 
line shows the cloud outline (cloud water fraction - 10 4 g g '). The surface gust front appears 
at —370 km. 


eration (Fovell et al. 1992) in this tropical case, perhaps 
due to the shallower penetration of the convection rel- 
ative to the higher tropical tropopause. 

As in the midlatitude simulations, the control case 


Surface Rain Rate 



Fig. 4. Contours of surface rainfall vs time and longitude (x) for 
the control simulation. The contour interval is 10 4 g g 1 (rainwater 
fraction). Intense rainfall is associated with the leading edge of the 
line (which moves westward) and rearward propagating cells (east- 
ward). 


(Fig. 3) has waves propagating both ahead of and to the 
rear of the westward moving squall line, but there is a 
clear dominance of rearward traveling gravity waves. 
The dominance of rearward propagating waves has re- 
cently been confirmed in observations of inertia-gravity 
waves from tropical Australia (Hamilton and Vincent 
1995). The easterly and westerly QBO phase cases (Fig. 
5) show dramatic evidence of wave absorption in shear 
zones and of ray path curvature owing to Doppler shift- 
ing of wave frequencies by mean flow shear. 

The damping of waves as they approach their critical 
levels has a profound influence on the wave spectra. 
Figure 6 shows the zonal phase speed spectra of mo- 
mentum flux for the westerly QBO case at several levels 
between 15 and 25 km. These spectra represent time 
averages over the 7-h analysis window. Near the tro- 
popause there is a broad distribution of waves in the 
phase speed range of ~±30 m s _1 . As height increases, 
there is an overall decrease in amplitude owing to prop- 
agation of waves out of the domain through the side 
boundaries. But between 20 and 24 km, in the strong 
westerly shear zone, there is a rapid filtering of the 
positive phase speed waves as waves of increasing phase 
speed encounter their critical levels. 

The effects of this filtering are clearly shown in the 
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Fig. 6. Spectra of momentum flux vs phase speed c in the - 16 m 
s ~ 1 moving frame. This is the intrinsic frame of reference at the 
tropopause. Spectra at z - 16, 18, 20, and 24 km are shown for the 
westerly shear case. Fluxes are computed as density times the cos- 
pectrum of u and vv' per unit phase speed. These spectra would he 
constant in height for linear, conservative waves in a closed domain, 
so changes reflect dissipation and cross-boundary fluxes. 


zonal- and time-averaged vertical profiles of eastward 
and westward momentum fluxes shown in Fig. 7. In the 
control case, the eastward momentum flux dominates 
over westward flux, and both the eastward and westward 
momentum fluxes decrease with height fairly uniformly. 
As will be shown in the next section, below the sponge 
layer this decrease with height is accounted for by waves 
leaving the domain at the lateral boundaries. In the east- 
erly QBO phase case the eastward momentum decreases 
substantially in the 15-20 km layer owing to the west- 
erly shear associated with the w-esterly wind peak at 20 
km. But between 21 and 24 km there is a very rapid 
decrease of westward momentum fluxes, corresponding 
to filtering by the main easterly shear zone. The situation 
is basically reversed in the westerly QBO phase case. 
In this case the westward momentum fluxes decrease 
substantially in the weak easterly shear zone below 20 
km, and the eastward propagating waves are nearly all 
removed between 20 and 23 km in the main westerly 
shear zone. 

4. Mean flow forcing 

The mean-flow forcing cannot be estimated simply 
by evaluating the net momentum flux convergence in 
the model domain because, as is clear from Figs. 3 and 

5, westward and eastward propagating waves both prop- 
agate horizontally out of the domain. The vertical mo- 
mentum fluxes computed over the model domain shown 
in Fig. 7 must be corrected for this effect. Such a cor- 
rection can be carried out using the wave activity con- 
servation principle. 

The method can be illustrated by considering a mon- 
ochromatic wave of phase 4> = Lx + mz — (ot, where 
k and m are the zonal and vertical wavenumbers, and 
ay is the frequency. Then for a slowly varying basic 


Control 
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Fig. 7. Three panels showing separate profiles of eastward and 
westward gravity wave momentum flux averaged over the (x,t) do- 
main for the control (top), easterly QBO (middle), and westerly QBO 
(bottom) cases. 


state, wave activity conservation for steady (i.e., time 
averaged) conservative waves can be expressed ap- 
proximately as (Andrews et al. 1987) 


d < c y » 4 ) + c>{c k: a) = o 
civ dz 


0) 


where angle brackets denote a time average, A = kE/ 
(oj — kU) is the “wave activity”(sometimes referred to 
as pseudomomentum), and the wave energy E is defined 
by 


E = 


A 


+ (gfl'/fl ,) 2 
2 2 2N 1 


( 2 ) 
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In (2), m\ vv\ and 6' denote the disturbance zonal ve- 
locity, vertical velocity, and potential temperature, re- 
spectively; overbars denote averages over a wave cycle; 
p s (z) and 0 s {z) denote basic state density and potential 
temperature, respectively; and N is the buoyancy fre- 
quency. Since, as shown below, c HZ A is equivalent to the 
vertical momentum flux, (1) can be used to estimate the 
momentum flux convergence associated with wave ac- 
tivity leaving the domain horizontally. 

Equation (1) can be expressed in a convenient form 
for numerical evaluation by using the relationships 
among the wave fields obtained from linear wave theory. 
If wave solutions are specified of the form 

(u\ vv\ 6 f ) = Re[(w, h\ 6) exp (/<£)], (3) 

where h, vt\ and 0 denote the amplitudes of the pertur- 
bations, then it can easily be shown (e.g., Holton 1992) 
that 


u = — {mfk)w 


(4) 


and 

(g6/6 s ) = -[iN 2 /(w - k(J)]w , (5) 

where k, m , and co are related by the dispersion equation 
(o> — kU) 2 /N 2 — k 2 /(k 2 4- m 2 ). (6) 


The zonal and vertical components of the group velocity 
are thus 


/da* d(o\ (m> —k)Nm 
C ‘" tg: ~ ( dk ' dm) ~ (k 2 + m'Y 12 ' (?) 

From (4)— <6) the wave energy of (2) can be expressed as 

u 2 (k 1 4 m 2 


£ = Ps 


( 8 ) 


and using the definition of A the vertical momentum 
flux can then be expressed as 

P ' 7 * 7 = = (9) 

Similarly, it can be shown that 



Substituting from (9) and (10) into (1) and integrating 
over the horizontal domain 0 < < L x then yields the 

desired result 


~~(p s u'w') 

pjz 


’<£> 


<*> 


2 

»mi 

2 



Thus the convergence of momentum flux owing to wave 
activity propagating sideways out of the domain can be 
easily measured as the difference in the zonal kinetic 
energy of the waves at the left and right boundaries of 
the domain. 

Although (11) has been derived for a monochromatic 


Boundary Flux Effects 



(m/s/day) 

Fig. 8. The effects of wave energy fluxes leaving the side bound- 
aries on the apparent mean-flow forcing for the control case. The 
dotted line shows the force per unit mass estimated directly from the 
convergence of the momentum flux profiles in Fig. 7. The dashed 
line shows the contribution due to side boundary fluxes from (1 1). 
The solid line then shows the zonal force corrected for waves leaving 
the domain. For this control case with no shear, the solid line gives 
the magnitude of unresolved errors in these calculations. 


wave, it is simple to extend the result to a spectrum of 
waves simply by summing the amplitudes on the right- 
hand side over all the spectral components. Figure 8 
shows the result of such a calculation for the control 
case. The boundary effect produces an apparent east- 
ward forcing in the domain in this case. This is con- 
sistent with the fact that eastward propagating waves 
are more strongly excited than westward propagating 
waves in this case and the eastward and westward modes 
have similar magnitudes of zonal group velocity. When 
the boundary effects calculated from the right-hand side 
of (11) are removed from the momentum flux diver- 
gence, the net zonal forcing computed for the control 
case is less than 0.5 m s _l day \ except in the sponge 
layer near the top where rapid damping of the westward 
propagating waves (Fig. 7) causes a net westward force 
of nearly 0.7 m s 1 day 1 . The smallness of the zonal 
force in the control case is evidence that, as expected 
for this case, gravity waves propagate through the strat- 
ospheric portion of the model domain with very little 
dissipation or mean-flow interaction. 

Similar zonal force calculations for the QBO easterly 
phase and QBO westerly phase cases (Figs. 9 and 10), 
on the other hand, show much different profiles of zonal 
forcing. In these cases the apparent forcing from waves 
propagating out of the domain is generally small com- 
pared to the forcing due to wave absorption in the mean- 
wind shear zones. As expected, a strong eastward forc- 
ing is centered in the main shear zone of the westerly 
QBO phase case (Fig. 10), and a strong westward forc- 
ing is centered in the main shear zone of the easterly 
QBO phase case (Fig. 9). 

The distributions of zonal forcing shown in Figs. 9 
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Boundary Flux Effects 



Fig. 9. Same as Fig. 8 but for the easterly shear case. Here the 
side boundary flux effects (dashed line) are similar to the control 
case, but small compared to the forcing estimated from the momentum 
flux convergence, which shows a strong westward peak in the easterly 
shear at ~22 km. 


Boundary Flux Effects 



(m/s/day) 

Fig. 10. Same as Fig. 8 but for the westerly shear case. As in Fig. 
9, the apparent forcing associated with side boundary fluxes (dashed 
line) is still small but enhanced slightly in comparison to the control 
(Fig. 8) because of the ray bending effects of the shear. The correction 
reduces the magnitude of the peak eastward forcing in the westerly 
shear zone in this case by ~ 1 m s 1 day l . 


and 10 are reminiscent of the strong forcing in shear 
zones characteristic of the LH68 model of the QBO. It 
is thus of some interest to compare the actual zonal 
forcing profile with that derived from the LH68 critical 
level absorption parameterization. According to this pa- 
rameterization, the zonal forcing may be expressed as 

--Y(PM^) = +k(i/) exp[(z - Zo)/H\\^, (12) 

p, dz dz 

where B 0 (U) — B 0 (c) at z 0 mapped onto the background 
wind field U(z) and B 0 (c) is defined so that poB 0 (c)dc is 
the momentum flux carried by waves in the phase speed 
range c to c+dc evaluated at z 0 , where p 0 = p/z 0 ). The 
mapping of B 0 (c) onto the background wind is subject 
to the condition that U(z) must be outside the range of 
background wind speeds between z () and z — dz (see 
Dunkerton 1996). 

The flux spectra B 0 (c) at z 0 = 15.875 km for the 
westerly and easterly QBO phase cases are shown in 
Fig. 11. As expected, these are nearly identical since 
background flow conditions in the two simulations are 
identical below 16 km. The zonal force profiles resulting 
from mapping these spectra onto the background winds 
are compared with the actual zonal force profiles derived 
from the model fields in Fig. 12. The agreement in mag- 
nitude in the main shear zone above 20 km is good in 
both cases, but the parameterized forcing profile com- 
puted from (12) is shifted upward about 1-1.5 km from 
the actual simulated forcing profile, and the magnitude 
is somewhat larger, consistent with the inverse depen- 
dence on density in (12). 

It is not surprising that wave damping occurs in the 
model as waves approach their critical levels. The dis- 
persion relation in (6) indicates that vertical wavelength 
is proportional to the Doppler-shifted phase speed. As 
the magnitude of the Doppler-shifted phase speed de- 


creases, damping by fourth-order diffusion becomes in- 
creasingly important. This damping can be represented 
as a scale-dependent Rayleigh friction and Newtonian 
cooling with rate coefficient a = Km 4 , where in the 
present case K = 1.628 X 10 6 m 4 s~‘. The decay of 
momentum flux with height owing to this linear damp- 
ing can be represented following HL72 (see also Dunk- 
erton 1996) as 


— ru: 


PoB 0 (c, k) 


■exp 


I D(c\ k, z')dz' 


\dc\dk . (13) 


Here 
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Fig. 1 1. Zonal intrinsic phase speed spectra of momentum flux at 
z = 15.875 km for the westerly (solid) and easterly (dotted) QBO 
cases. 
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Fig. 12. Zonal force per unit mass estimated from the spectra in Fig. 1 1 and the LH68 model [Eq. (12)] plotted as dashed lines for (a) 
westerly and (b) easterly QBO cases. The boundary- flux-corrected forcing profiles from Figs. 10 and 9 are also plotted with solid lines for 
comparison. 


2 Km A N _ 2KN 5 
k(c - U) 1 ~~ k{c - t/) 6 ’ 


( 14 ) 


where the far right expression holds for m 2 k 2 . It is 
easy to see that in this case the damping for each wave 
component will be strongly peaked just below the crit- 
ical level as c — U — » 0. As shown in the appendix, an 
asymptotic analysis indicates that the maximum wave- 
induced force will then occur about 1.46 km below the 
critical level in the easterly shear case and 0.75 km 
below the critical level in the westerly shear case. The 
momentum flux convergence obtained from finite dif- 
ferencing the expression given in (13), shown in Fig. 
13, provides a good approximation to the actual flux 
convergence of the simulated waves. 

Hyperdiffusion in the model is analogous to the dis- 
sipative processes acting on gravity waves approaching 
critical levels in the real atmosphere. Convective and 
dynamical instability, and scale-selective radiative 


damping, are all likely to become important as the ver- 
tical wavelength decreases. Thus, as in the simulations, 
the zonal force produced by a given wave would be 
concentrated a small distance below the critical level. 
Linear saturation theory applied to waves with small 
amplitudes at their source levels would produce a similar 
result (e.g., Alexander and Rosenlof 1996). 

The mean zonal wind acceleration, evaluated by tak- 
ing the difference between a slightly smoothed version 
of the basic-state profile and the zonally averaged wind 
profile at hour 1 1 of the simulation divided by the time 
interval, is compared with the time-averaged wave-in- 
duced zonal force in Fig. 13. Although the mean wind 
tendencies in the main shear zones are in the same sense 
as the forcing for both the westerly and easterly QBO 
cases, the amplitude of the wind acceleration is much 
smaller than that of the wave-induced forcing. This fail- 
ure of wave-induced forcing to be completely balanced 
by local flow deceleration is not surprising for a model 
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Fig. 13. Zonal force per unit mass estimated from tropopause momentum flux spectra and the dissipation model (Eq. (13); dashed lines j 
for the (a) westerly and (b) easterly QBO cases. These are compared to the boundary-flux-corrected profiles from Figs. 10 and 9 (solid 
lines). Also plotted, as dotted lines, is the change in the mean zonal wind observed in the simulations between the initial time and t - 1 1 
h divided by the 1 1-h time interval. 
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with open lateral boundaries. Even though we have 
corrected the zonal force for the effects of gravity 
waves leaving the domain, there is still a net dynamic 
pressure gradient across the domain, which provides 
part of the momentum balance. As discussed by Durran 
(1995) in the context of mountain waves, the local 
wave-induced forcing generates pressure disturbances 
that propagate far from the forcing region and create 
a highly nonlocal mean flow response. The wave-in- 
duced force is thus not balanced by local accelerations, 
but only by the zonal flow acceleration averaged 
around a complete latitude circle. Hence, for evaluation 
of the possible impact of convectively generated grav- 
ity waves on the QBO, it is appropriate to use the 
computed wave-induced force, not the mean zonal 
wind change in the model domain. 

5. Discussion: Implications for the QBO 

In the observed equatorial QBO the descent rate of 
the westerly shear zone tends to be faster and more 
regular than the descent rate of the easterly shear zone. 
These differences are partly due to modulation of the 
mean upwelling rate by a secondary meridional circu- 
lation associated with the QBO and partly due to dif- 
fering distributions of eastward and westward wave- 
induced forcing. Neglecting the asymmetries between 
the two phases, the average zonal acceleration rate in 
the QBO shear zones is of order 6ms 1 month \ or 
0.2 m s 1 day 1 (Dunkerton and Delisi 1985), averaged 
over the area equatorward of about ±10° latitude. But, 
as noted by Dunkerton (1996), since the downward 
propagation of the QBO is opposed by the upwelling 
Brewer-Dobson circulation, a zonal force nearly twice 
the magnitude of the zonal wind acceleration is required 
to produce the observed oscillation. 

The gravity wave-induced zonal forces in the west- 
erly and easterly QBO phase squall line simulations 
shown in Figs. 9 and 10 are of order 5 m s _l day 1 
averaged over the model domain. If it is assumed that 
the length scale parallel to the simulated squall line is 
about 200 km, then the squall line domain is about 0.2% 
of the total area equatorward of ±10° latitude. Thus, if 
10 such squall lines were active on average over the 
entire tropical region, they could account for about 1/4 
of the required zonal forcing of the QBO. 

Climatological information on the areal extent of 
tropical convective systems has been reported by Zhang 
(1993) and Chen et al. (1996). According to S. Chen 
(1996, personal communication), the analysis of the lat- 
ter paper implies that 6% of the area of the Indian Ocean 
and western Pacific was covered by convective systems 
during the intensive observation period (November 
1992-February 1993) of the Tropical Ocean Global At- 
mosphere Coupled Ocean-Atmosphere Response Ex- 
periment (TOGA COARE). Since TOGA COARE ob- 
servations were for a region and time of very active 


convection, the above estimate of 2% coverage for the 
entire equatorial area appears to be reasonable. 

The West African squall line simulated in this study 
occurs, however, in an environment of unusually high 
conditional instability, so it is likely to have more in- 
tense convective cells than typical of convective storms 
over most of the Tropics. Also, squall lines in nature 
are generally not aligned strictly in the east-west di- 
rection, but at some angle, and this would further reduce 
the estimated zonal momentum carried by the gravity 
waves. Furthermore, preliminary studies with a three- 
dimensional version of the cloud model suggest that 
simulated gravity wave amplitudes are larger in the two- 
dimensional simulation than in three-dimensions. Thus, 
it is likely that mesoscale gravity waves generated by 
such convection actually account for somewhat less than 
1/4 of the acceleration in the QBO. 

Nevertheless, mesoscale gravity waves may play an 
important role in accounting for the observed lateral 
width of the equatorial QBO. LH68 focused on large- 
scale equatorially trapped wave modes, rather than me- 
soscale gravity waves, partly because such modes are 
trapped in the equatorial region and hence can appar- 
ently readily account for the observed decay of the QBO 
outside of the equatorial zone. Previous two-dimen- 
sional and three-dimensional simulations of the QBO, 
however, have tended to produce oscillations with me- 
ridional scales smaller than observed (e.g.. Plumb and 
Bell 1982). Any model that depends solely on equatorial 
wave modes is likely to have this deficiency, since the 
meridional width of an equatorial mode shrinks rapidly 
as the wave approaches its critical level. Thus, what 
must be explained is not why the QBO wind oscillation 
is confined to the equatorial zone, but rather why the 
QBO is not more tightly confined to the equatorial zone 
than is observed. 

Since vertically propagating mesoscale gravity 
waves will be generated wherever convection occurs, 
such waves should be distributed over a wide range of 
latitudes, rather than being confined to the equatorial 
zone. If such waves provide a substantial fraction of 
the QBO zonal forcing, it might be argued that the 
QBO should not be equatorially confined. P. Haynes 
(1994, personal communication) has shown, however, 
that the response to a wave-induced zonal force in the 
stratosphere depends strongly on the latitude where the 
force is applied. At low latitudes the force is balanced 
by a zonal wind acceleration, while at higher latitudes 
it is balanced by the Coriolis force of an induced me- 
ridional circulation [the extratropical pumping dis- 
cussed in Holton et al. (1995)]. Consequently, even if 
the wave-induced forcing were distributed uniformly 
in latitude, the resulting zonal acceleration should be 
equatorially trapped. Thus, it is tempting to conclude 
that mesoscale gravity waves do have a significant role 
to play in the equatorial QBO, especially in accounting 
for the observed wind oscillation several degrees of 
latitude away from the equator. 
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zonal wavelength of 62.8 km the maximum wave-in- 
duced force thus occurs at ( z c — z) = 146 km in the 
easterly QBO case and ( z c — z) — 0.75 km in the west- 
erly QBO case. 
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APPENDIX 

Zonal Force Induced by Diffusion 

The distance below the critical level at which the 
maximum wave-induced force occurs for waves subject 
to the fourth-order vertical diffusion in the numerical 
model can be computed using the expressions derived 
in the appendix of Dunkerton (1996). From (13) the 
momentum flux per unit phase speed per unit wave- 
number is 


p 0 B(c, k \ z) = poB 0 (c , k) exp 
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D(c, k , z ') dz 


, (Al) 


with D(c, k, z) given by (14). Following Dunkerton the 
argument in the exponent in (Al) can be expressed as 
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where A = {dU/dz) Zc . Because the basic-state shear is 
constant in the main shear layers for both the easterly 
shear and westerly shear simulations, c - U = A (z c ~ 
z). Thus, (A2) can be expressed as 


£ 


D dz 


2KAY 


5kA 6 (z r - z ) 5 


€ 5 


(A3) 


and D can be expressed as 
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where £ = e(z c - zY ' and 6 = N[2K/(5 k A 6 )] 1 ' 5 . 

The wave-induced force per unit mass per unit phase 
speed per unit wavenumber is 
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I d Po 

---[ Po B(c, k, z)] — —DB 
Pv dz A 


= exp(-f 5 ) expf^*). (A5) 

This has its maximum for £ — 1.04. In the simulations 
reported here K = 1 .628 X 10 6 m 4 s' 5 and N = 2.2X 10~ 2 
s _i . For the westerly shear case A = 7. OX 10" 3 m _l and 
for the easterly shear case A = -4. OX 10 -3 m _l . For a 
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